*** One-Child Policy, Marriage Distortion, and Welfare Loss ***
*** by Wei Huang, Yinghao Pan and Yi Zhou ***

clear all 
set more off
set scheme s1mono
cap log close 
graph set window fontface "Times New Roman"

set niceness 10
gl MAIN = "F:\Dropbox"

gl path1 = "$MAIN/Census"
gl path2 = "$MAIN/Marriage/Working_data"
gl path3 = "$MAIN/Marriage/Restat/dataverse_files/Results" 
gl path4 = "$MAIN/Marriage/Restat/dataverse_files/Graphics" 




**** Unconditional rates ****

use "$path2/unmarried", clear
keep if relation == 0 | relation == 1 | relation == 2 | relation == 6 

merge 1:1 wave area hhid name relation using "$path2/min_marr_policy", keepusing(marry_1 marry_2 marry_3) nogen

forvalues i = 30(1)35{
gen year = year_birth + `i' 
merge m:1 prov year using "$path2/fines_use", keepusing(fine)

drop if _merge == 2
replace fine = 0 if fine == . & year <= 1979
drop _merge year
ren fine fine_`i'
}

drop if mi(hktype)
replace urban = hktype == 1 
su age 

cap gen female = 1-male
cap gen han = 1- minority 
cap gen se = .
replace unmarried = unmarried * 100 
gen unmarr_25 = marr_age > 25 if !mi(marr_age)
replace unmarr_25 = 1 if unmarried == 100 
replace unmarr_25 = unmarr_25 *100 

gen urban_hk = hktype == 2 


gen marry_9 = unmarried == 0 & (mi(marry_1) & mi(marry_2) & mi( marry_3))

recode marry_1 marry_2 marry_3 (. = 0) 

gen HH = marry_1 * 100 
gen HM = marry_2 * 100 
gen MM = marry_3 * 100 
gen Other = marry_9 * 100 
gen HOMO = HH + MM
su HOMO

keep if age <= 55 
keep if age > 25
keep if year_birth <= 1980 
keep if Other == 0 



****  Figure 1(a) ****

preserve
keep if wave >= 2005 
keep if han == 1 
cap gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | prov == 35 ///
|prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
gen category = 1 if han == 0 & treat == 0 
replace category = 2 if han == 0 & treat == 1 

gen treatXyear = treat* year_birth 
recode treatXyear (1965 = 0), gen(treatXDyear)


recode year_birth (1900/1965 = 1) (1966/1970 = 2) (1971/1980 = 3), gen(year_birth_cate)
gen treatXcate = treat* year_birth_cate 

reghdfe unmarried ib1.treatXcate , a(year_birth treat)   cl(year_birth)
test 3.treatXcate
loc F: di %5.3g r(F)
loc p: di %3.2f r(p)
restore

preserve 
keep if wave >= 2005
keep if han == 1 
cap gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | prov == 35 ///
|prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
gen category = 1 if han == 1 & treat == 0 
replace category = 2 if han == 1 & treat == 1 
replace category = 3 if han == 0 
gen treatXyear = treat* year_birth 


levelsof year_birth
qui:{
foreach x of numlist `r(levels)' {

ttest unmarried if year_birth==`x' ,by(category)
replace se=r(se) if year_birth==`x'
}
}



collapse unmarried se, by(year_birth category)

reshape wide unmarried , i(year_birth) j(category)

foreach var in "unmarried"{
gen `var'_diff = `var'2 - `var'1
gen ub = `var'_diff + 1.96* se
gen lb = `var'_diff - 1.96* se

}
reg unmarried_diff year_birth if year_birth < 1965, cl(year_birth)
testparm year_birth
loc F0: di %5.3g r(F)
loc p0: di %3.2f r(p)

tw (con unmarried1 unmarried2 year if year_birth > 1940, xlabel(1950(5)1980) xtit("Birth cohort") ytit("Unmarried (%)") m(O D)) (line unmarried_diff lb ub year, yaxis(2) ytit("Difference (%)", axis(2))  lp(dash dot dot) lw(thick ) lc(gray gray gray) ///
legend(order(2 1 3) lab(1 "Provinces w/o fines increase in 1989-1995") lab(2 "Provinces w/ fines increase in 1989-1995") lab(3 "Difference (%)") ///
size(small) pos(11) ring(0) col(1)) ylabel(0(3)24, grid) ylabel(-2(1)6, axis(2)) m(O D) yline(0, axis(2) lp(dash) lc(gray))  xline(1965 1970, lp(dash)) ///
text(18 1965 "Aged 25 ""in 1990", size(small)) text(18 1970 "Aged 20 ""in 1990", size(small)) ///
text(12 1955 "Pretrend test" "F-stat =`F0'" "p-val = `p0'", size(small)) text(2 1975 "DID test""F-stat = `F'" "p-value = `p'", size(small))) 
gr export "$path4/unmarried_han_cohort.pdf",replace 
restore 


****  Figure 1(b) ****

preserve
keep if wave >= 2005 
keep if han == 0 
cap gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | prov == 35 ///
|prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
gen category = 1 if han == 0 & treat == 0 
replace category = 2 if han == 0 & treat == 1 
gen treatXyear = treat* year_birth 




recode year_birth (1900/1965 = 1) (1966/1970 = 2) (1971/1980 = 3), gen(year_birth_cate)

gen treatXcate = treat* year_birth_cate 

reghdfe unmarried ib1.treatXcate , a(year_birth treat) cl(year_birth)
test 3.treatXcate
loc F: di %3.2f r(F)
loc p: di %3.2f r(p)
restore

preserve 
keep if wave >= 2005 
keep if han == 0 
cap gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | prov == 35 ///
|prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
gen category = 1 if han == 0 & treat == 0 
replace category = 2 if han == 0 & treat == 1 
drop if mi(category)

levelsof year_birth
foreach x of numlist `r(levels)' {
ttest unmarried if year_birth==`x' ,by(category)
replace se=r(se) if year_birth==`x'
}


collapse unmarried se, by(year_birth category)

reshape wide unmarried , i(year_birth) j(category)


foreach var in "unmarried"{
gen `var'_diff = `var'2 - `var'1
gen ub = `var'_diff + 1.96* se
gen lb = `var'_diff - 1.96* se
}

reg unmarried_diff year_birth if year_birth < 1965, cl(year_birth)
testparm year_birth
loc F0: di %5.3g r(F)
loc p0: di %3.2f r(p)

tw (con unmarried1 unmarried2 year if year_birth > 1940, xlabel(1950(5)1980) xtit("Birth cohort") ytit("Unmarried (%)") m(O D)) (line unmarried_diff lb ub year, yaxis(2) ytit("Difference (%)", axis(2))  lp(dash dot dot) lw(thick ) lc(gray gray gray) ///
legend(order(2 1 3) lab(1 "Provinces w/o fines increase in 1989-1995") lab(2 "Provinces w/ fines increase in 1989-1995") lab(3 "Difference (%)") ///
size(small) pos(11) ring(0) col(1)) ylabel(0(3)24, grid) ylabel(-4(1)6, axis(2)) yline(0, axis(2) lp(dash) lc(gray)) m(O D) xline(1965 1970, lp(dash)) ///
text(18 1965 "Aged 25 ""in 1990", size(small)) text(18 1970 "Aged 20 ""in 1990", size(small)) ///
text(15 1955 "Pretrend test"  "F-stat =`F0'" "p-val = `p0'", size(small)) text(2 1975 "DID test""F-stat = `F'" "p-val = `p'", size(small)))
gr export "$path4/unmarried_mino_cohort.pdf",replace 

restore 

****  Figure 2(a) ****


preserve
keep if unmarried == 0 
keep if wave >= 2005 
keep if mino_policy == 1
cap gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | prov == 35 ///
|prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
gen category = 1 if han == 0 & treat == 0 
replace category = 2 if han == 0 & treat == 1 
gen treatXyear = treat* year_birth 
recode year_birth (1900/1965 = 1) (1966/1970 = 2) (1971/1980 = 3), gen(year_birth_cate)

gen treatXcate = treat* year_birth_cate 

reghdfe HOMO ib1.treatXcate , a(year_birth treat) cl(year_birth) 
test 3.treatXcate
loc F: di %5.3g r(F)
loc p: di %3.2f r(p)
restore


preserve
keep if unmarried == 0 
keep if wave >= 2005
keep if han
keep if mino_policy == 1
cap gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | prov == 35 ///
|prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
gen category = 1 if treat == 0 & mino_policy == 1 
replace category = 2 if  treat == 1 & mino_policy == 1 


drop if mi(category)

levelsof year_birth
foreach x of numlist `r(levels)' {
ttest HOMO if year_birth==`x' ,by(category)
replace se=r(se) if year_birth==`x'

}

collapse HOMO se, by(year_birth category)

recode year_birth (1930/1969 = 1) (1970/1990 = 2), gen(yob_group)


reshape wide HOMO , i(year_birth) j(category)

foreach var in "HOMO"{
gen `var'_diff = `var'2 - `var'1
gen ub = `var'_diff + 1.96* se
gen lb = `var'_diff - 1.96* se
}

reg HOMO_diff year_birth if year_birth < 1965, cl(year_birth)
testparm year_birth
loc F0: di %5.3g r(F)
loc p0: di %3.2f r(p)

tw (con HOMO1 HOMO2 year_birth if year_birth > 1940, xlabel(1950(5)1980) xtit("Birth cohort") ytit("Intra-ethnic marriage (%), conditional") m(O D)) ///
 (line HOMO_diff  lb ub year, yaxis(2) ytit("Difference (%)", axis(2))  lp(dash dot dot) lw(thick ) lc(gray gray gray) ///
legend(order(2 1 3) lab(1 "Provinces w/o fines increase in 1989-1995") lab(2 "Provinces w/ fines increase in 1989-1995") lab(3 "Difference (%)") ///
size(small) pos(7) ring(0) col(1)) ylabel(92(1)100, grid) ylabel(-4(1)4, axis(2)) m(O D) yline(0, axis(2) lp(dash) lc(gray)) ///
 xline( 1965.5 1970, lp(dash))  text(99 1965 "Aged 25 ""in 1990", size(small)) text(99 1970 "Aged 20 ""in 1990", size(small))  ///
text(94.2 1953 "Pretrend test"  "F-stat =`F0'" "p-val = `p0'", size(small)) text(93 1975 "DID test""F-stat = `F'" "p-value = `p'", size(small)))

 gr export "$path4/homo_policy_cohort.pdf",replace 
restore
 
 

****  Figure 2(b) ****


preserve
keep if unmarried == 0 
keep if wave >= 2005 
keep if mino_policy == 0
cap gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | prov == 35 ///
|prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
gen category = 1 if han == 0 & treat == 0 
replace category = 2 if han == 0 & treat == 1 
gen treatXyear = treat* year_birth 
recode year_birth (1900/1965 = 1) (1966/1970 = 2) (1971/1980 = 3), gen(year_birth_cate)

gen treatXcate = treat* year_birth_cate 

reghdfe HOMO ib1.treatXcate , a(year_birth treat) cl(year_birth) 
test 3.treatXcate

loc F_no: di %3.2f r(F)
loc p_no: di %3.2f r(p)
restore



 preserve
keep if unmarried == 0 
keep if wave >= 2005
keep if han
keep if mino_policy == 0
cap gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | prov == 35 ///
|prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
gen category = 1 if treat == 0  
replace category = 2 if  treat == 1 

levelsof year_birth
foreach x of numlist `r(levels)' {
ttest HOMO if year_birth==`x' ,by(category)
replace se=r(se) if year_birth==`x'

}

collapse HOMO se, by(year_birth category)


reshape wide HOMO , i(year_birth) j(category)

foreach var in "HOMO"{
gen `var'_diff = `var'2 - `var'1
gen ub = `var'_diff + 1.96* se
gen lb = `var'_diff - 1.96* se
}
reg HOMO_diff year_birth if year_birth < 1965, cl(year_birth)
testparm year_birth
loc F0: di %5.3g r(F)
loc p0: di %3.2f r(p)

 tw (con HOMO1 HOMO2 year if year_birth > 1940, xlabel(1950(5)1980) xtit("Birth cohort") ytit("Intra-ethnic marriage (%), conditional") m(O D)) ///
 (line HOMO_diff  lb ub year, yaxis(2) ytit("Difference (%)", axis(2))  lp(dash dot dot) lw(thick ) lc(gray gray gray) ///
legend(order(2 1 3) lab(1 "Provinces w/o fines increase in 1989-1995") lab(2 "Provinces w/ fines increase in 1989-1995") lab(3 "Difference (%)") ///
size(small) pos(7) ring(0) col(1)) ylabel(92(1)100, grid) ylabel(-4(1)4, axis(2)) m(O D)  yline(0, axis(2) lp(dash) lc(gray)) ///
 xline( 1965 1970, lp(dash))  text(98 1965 "Aged 25 ""in 1990", size(small)) text(98 1970 "Aged 20 ""in 1990", size(small)) ///
 text(95 1955 "Pre-trend test""F-stat = `F0'" "p-value = `p0'", size(small)) text(93 1975 "DID test""F-stat = `F_no'" "p-value = `p_no'", size(small)))

gr export "$path4/homo_nonpolicy_cohort.pdf",replace 

restore 




********* Appendix Figure A4 ***********

**** Appendix Figure A4(a) ****
set scheme s1color

preserve
keep if wave >= 2005 

recode year_birth (1900/1965 = 1) (1966/1970 = 2) (1971/1980 = 3), gen(year_birth_cate)

gen treatXcate = mino_policy* year_birth_cate 

reghdfe HOMO ib1.treatXcate , a(year_birth mino_policy) cl(year_birth) 
test 3.treatXcate
loc F_uncon: di %5.3g r(F)
loc p_uncon: di %3.2f r(p)
restore

preserve

keep if wave >= 2005

levelsof year_birth
foreach x of numlist `r(levels)' {
ttest HOMO if year_birth==`x' ,by(mino_policy)
replace se=r(se) if year_birth==`x'

}

collapse HOMO fine_16_25 se, by(mino_policy year_birth)
reshape wide HOMO fine_16_25, i(year_birth) j(mino_policy)
foreach var in "HOMO"{
gen `var'_diff = `var'1 - `var'0
gen ub = `var'_diff + 1.96* se
gen lb = `var'_diff - 1.96* se
}

reg HOMO_diff year_birth if year_birth < 1965, cl(year_birth)
testparm year_birth
loc F0: di %5.3g r(F)
loc p0: di %3.2f r(p)

tw (con HOMO0 HOMO1 year_birth,  ylabel(80(4)100) xtit(Birth cohort) ytit("Intra-ethnic marriage (%, unconditional)") m(O T) ylabel(72(4)100, grid)) ///
(line HOMO_diff  lb ub year, yaxis(2) ytit("Difference (%)", axis(2))  lp(dash dot dot) lw(thick ) lc(blue blue blue) ///
 ylabel(-8(2)-2, axis(2)) ///
legend(order(1 "No Pref-policy regions" 2 "Pref-policy regions" 3 "Difference (%)") size(small)  ring(0) col(1) pos(7)) ///
xline( 1965 1970, lp(dash)) text(85 1965 "Aged 25 ""in 1990") text(85 1970 "Aged 20 ""in 1990")  text(75 1974 "F-stat = `F_uncon'" "p-value = `p_uncon'"))
gr export "$path4/homo_policy.pdf", replace 

restore 



**** Appendix Figure A4(b) ****



preserve
keep if unmarried == 0
keep if wave >=2005 

recode year_birth (1900/1965 = 1) (1966/1970 = 2) (1971/1980 = 3), gen(year_birth_cate)

gen treatXcate = mino_policy* year_birth_cate 

reghdfe HOMO ib1.treatXcate , a(year_birth mino_policy) cl(year_birth) 
test 3.treatXcate

loc F_con: di %5.3g r(F)
loc p_con: di %3.2f r(p)
restore



preserve
keep if unmarried == 0
keep if wave >=2005 

levelsof year_birth
foreach x of numlist `r(levels)' {
ttest HOMO if year_birth==`x' ,by(mino_policy)
replace se=r(se) if year_birth==`x'

}

collapse HOMO fine_16_25 se, by(mino_policy year_birth)
reshape wide HOMO fine_16_25, i(year_birth) j(mino_policy)
foreach var in "HOMO"{
gen `var'_diff = `var'1 - `var'0
gen ub = `var'_diff + 1.96* se
gen lb = `var'_diff - 1.96* se
}
tw (con HOMO0 HOMO1 year_birth,  ylabel(90(2)100) xtit(Birth cohort) ytit("Intra-ethnic marriage (%, conditional)") m(O T) ylabel(90(2)100, grid)) ///
(line HOMO_diff   lb ub year, yaxis(2) ytit("Difference (%)", axis(2))  lp(dash dot dot) lw(thick ) lc(blue blue blue) ///
 ylabel(-8(2)-2, axis(2)) ///
legend(order(1 "No Pref-policy regions" 2 "Pref-policy regions" 3 "Difference (%)") size(small)  ring(0) col(1) pos(7)) ///
xline( 1965 1970, lp(dash)) text(93 1965 "Aged 25 ""in 1990") text(93 1970 "Aged 20 ""in 1990")  text(91 1975 "F-stat = `F_con'" "p-value = `p_con'"))
gr export "$path4/homo_policy_condition.pdf", replace 

restore 



set scheme s1mono




keep if sex == 2 // women 

egen n_birth = rowtotal(numbirm numbirf)
replace n_birth = . if age > 45

****  Figure 3(a) ****

preserve

keep if wave >= 2005 
keep if han == 1 
cap gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | prov == 35 ///
|prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
gen category = 1 if han == 0 & treat == 0 
replace category = 2 if han == 0 & treat == 1 
gen treatXyear = treat* year_birth 
recode year_birth (1900/1965 = 1) (1966/1970 = 2) (1971/1980 = 3), gen(year_birth_cate)

gen treatXcate = treat* year_birth_cate 

reghdfe n_birth ib1.treatXcate , a(year_birth treat) cl(year_birth) 
test 3.treatXcate 

loc F_han: di %5.3g r(F)
loc p_han: di %3.2f r(p)
restore




preserve 
keep if wave == 2005
keep if han == 1 
gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | ///
prov == 35 |prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
keep if !mi(n_birth)
keep if year_birth >= 1960

gen category = 1 if han == 1 & treat == 0 
replace category = 2 if han == 1 & treat == 1 
replace category = 3 if han == 0 


levelsof year_birth
foreach x of numlist `r(levels)' {
ttest n_birth if year_birth==`x' ,by(category)
replace se=r(se) if year_birth==`x'

}


collapse n_birth se, by(year_birth category)

reshape wide n_birth , i(year_birth) j(category)

foreach var in "n_birth"{
gen `var'_diff = `var'2 - `var'1
gen ub = `var'_diff + 1.96* se
gen lb = `var'_diff - 1.96* se
}
reg  n_birth_diff year_birth if year_birth <= 1970 , cl(year_birth)
test year_birth
loc F0: di %5.3g r(F)
loc p0: di %3.2f r(p)

tw (con n_birth1 n_birth2 year if year_birth >= 1960, xtit("Birth cohort") xlabel(1960(5)1980) ytit("# of births") m(O D)) (line n_birth_diff lb ub year, yaxis(2) ytit("Difference", axis(2))  lp(dash dot dot) lw(thick ) lc(gray gray gray) ///
legend(order(2 1 3) lab(1 "Provinces w/o fines increase in 1989-1995") lab(2 "Provinces w/ fines increase in 1989-1995") lab(3 "Difference") ///
size(small) pos(1) ring(0) col(1))  ylabel(0(.5)3, grid) ylabel(-0.2(.1).3, axis(2))  m(O D) xline( 1965 1970, lp(dash)) yline(0, axis(2) lp(dash) lc(gray)) text(1 1965 "Aged 25 ""in 1990", size(small)) text(1 1970 "Aged 20 ""in 1990", size(small)) ///
 text(0.5 1962 "Pretrend test""F-stat = `F0'" "p-value = `p0'", size(small))  text(0.5 1975 "DID test" "F-stat = `F_han)'" "p-value = `p_han'", size(small)))
gr export "$path4/n_birth_cohort.pdf", replace 
restore 

****  Figure 3(b) ****
preserve

keep if wave >= 2005 
keep if han == 0 
cap gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | prov == 35 ///
|prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
gen category = 1 if han == 0 & treat == 0 
replace category = 2 if han == 0 & treat == 1 
gen treatXyear = treat* year_birth 
recode year_birth (1900/1965 = 1) (1966/1970 = 2) (1971/1980 = 3), gen(year_birth_cate)

gen treatXcate = treat* year_birth_cate 

reghdfe n_birth ib1.treatXcate , a(year_birth treat) cl(year_birth) 
test 3.treatXcate 

loc F_min: di %3.2f r(F)
loc p_min: di %3.2f r(p)
restore
preserve 
keep if wave == 2005
keep if han == 0 
gen treat = 1 if prov == 11 | prov == 13 | prov == 21 | prov == 31 | prov == 32 | prov == 33 | prov == 35 |prov == 41 | prov == 42 | prov == 43 |prov == 44 | prov == 45| prov == 46| prov == 53| prov == 64 | prov == 65
replace treat = 0 if mi(treat)
drop if treat == . 
keep if !mi(n_birth)
keep if year_birth >= 1960

gen category = 1 if han == 0 & treat == 0 
replace category = 2 if han == 0 & treat == 1 

levelsof year_birth
foreach x of numlist `r(levels)' {
ttest n_birth if year_birth==`x' ,by(category)
replace se=r(se) if year_birth==`x'

}


collapse n_birth se, by(year_birth category)


reshape wide n_birth , i(year_birth) j(category)

foreach var in "n_birth"{
gen `var'_diff = `var'2 - `var'1
gen ub = `var'_diff + 1.96* se
gen lb = `var'_diff - 1.96* se
}
reg  n_birth_diff year_birth if year_birth <= 1970 , cl(year_birth)
test year_birth
loc F0: di %5.3g r(F)
loc p0: di %3.2f r(p)


tw (con n_birth1 n_birth2 year if year_birth >= 1940, xtit("Birth cohort") ytit("# of births") m(O D)) (line n_birth_diff lb ub year, yaxis(2) ytit("Difference", axis(2))  lp(dash dot dot) lw(thick ) lc(gray gray gray) ///
legend(order(2 1 3) lab(1 "Provinces w/o fines increase in 1989-1995") lab(2 "Provinces w/ fines increase in 1989-1995") lab(3 "Difference") ///
size(small) pos(1) ring(0) col(1)) ylabel(0.5(.5)3, grid) ylabel(-0.4(.2).6, axis(2)) m(O D) xline( 1965 1970, lp(dash)) yline(0, axis(2) lp(dash) lc(gray)) ///
text(1 1965 "Aged 25 ""in 1990", size(small)) text(1 1970 "Aged 20 ""in 1990", size(small)) ///
 text(0.65 1962 "Pretrend test""F-stat = `F0'" "p-value = `p0'", size(small)) ///
 text(0.65 1978 "F-stat = `F_min)'" "p-value = `p_min'", size(small)))
gr export "$path4/n_birth_cohort_mino.pdf", replace 
restore 



****  Figure 5 ****



use "$path2/min_marr_policy", clear 
drop if mi(hktype)
su age 
egen sex_age = group(male age)
egen hukou_prov = group(prov_hk hktype)
egen mino_birth = group(minority year_birth)
gen urban_hk = hktype == 2 
cap gen female = 1-male
cap gen han = 1- minority 
replace children = . if age > 45 // NOTE THIS 

keep if age <= 45
keep if age > 25
keep if year_birth <= 1980


cap drop birth_g sex_birth_g group

egen group = group(hukou_prov)
egen awt = sum(aw), by(group han)
gen aw_m = aw*marry_2
egen awm = sum(aw_m), by(group han)
gen age_g = int(age/5)
gen birth_g = int(year_birth/10)
tab age_g
tab birth_g

replace age_g = 10 if age_g == 11

gl CONTROL_S = "  i.birth_g i.age_g#sex"
drop if mi(group )
gen coef_h_m = . 
gen coef_h_c = . 
gen coef_m_m = . 
gen coef_m_c = . 

egen fine_16_25 = rowmean(fine_16-fine_25)
qui:{
forvalues i = 1(1)67{
cap reg marry_2 fine_16_25 $CONTROL_S [aw = aw] if group == `i' & han == 1 
cap replace coef_h_m = _b[fine_16_25] if group == `i' & han == 1 
cap reg children fine_16_25 $CONTROL_S [aw = aw] if group == `i' & marry_2 == 1  & han == 1 
cap replace coef_h_c = _b[fine_16_25] if group == `i' & han == 1 
}
}

fcollapse  coef_* mino_policy awt awm, by(group prov_hk hktype han)
duplicates drop 

replace coef_h_m = coef_h_m * 100 

keep if han == 1 
	  
cap drop outlier
gen outlier = . 
replace outlier = 0 if outlier == . 

replace outlier = 1 if abs(coef_h_m) >= 2 
replace outlier = 1 if abs(coef_h_c) >= 2 

****  Figure 5(a) ****

													
reg coef_h_c coef_h_m [aw = awt] if  mino_policy <0.05 & !outlier, r
local beta = string(_b[coef_h_m],"%9.3f")
local cons = string(_b[_cons],"%9.3f")
local se = string(_se[coef_h_m], "%9.3f")
tw (scatter  coef_h_c coef_h_m  if  mino_policy <0.05 & !outlier [aw = awt], ms(Oh) xlabel(-1.5(0.5)1.5) xline(0, lp(dash)) yline(0, lp(dash)) ylabel(-2(0.5)1)) ///
(lfit coef_h_c coef_h_m if  mino_policy <0.05 & !outlier  [aw = awt] ,legend(off) xtit("Impact on Inter-ethnic Marriage") ytit("Impact on # of Births in Inter-ethnic Marriage") ///
text(-1.5 0.5 "y =  `beta'x `cons'" "(`se')" ))
gr export "$path4/fig_imp_marriage_child0.pdf",replace 
	
****  Figure 5(b) ****

reg coef_h_c coef_h_m [aw = awt] if  mino_policy >0.5  & !outlier, r
local beta = string(_b[coef_h_m],"%9.3f")
local cons = string(_b[_cons],"%9.3f")
local se = string(_se[coef_h_m], "%9.3f")
tw (scatter  coef_h_c coef_h_m if  mino_policy >0.5  & !outlier [aw = awt], msymbol(Oh)  xlabel(-1.5(0.5)1.5)  xline(0, lp(dash)) yline(0, lp(dash)) ylabel(-2(0.5)1)) ///
(lfit coef_h_c coef_h_m [aw = awt] if  mino_policy >0.5 & !outlier  , legend(off) xtit("Impact on Inter-ethnic Marriage") ytit("Impact on # of Births in Inter-ethnic Marriage") ///
text(-1.5 0.5 "y =  `beta'x `cons'" "(`se')") )
gr export "$path4/fig_imp_marriage_child1.pdf",replace 

gen min = mino_policy >0.5 & !outlier
reg coef_h_c c.coef_h_m##i.min [aw = awt] if  (mino_policy >0.5  & !outlier) | (mino_policy <0.05 & !outlier), r
test c.coef_h_m # 1.min = c.coef_h_m # 0.min

